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Abstract. 

A discrete asymmetric exclusion process (ASEP) is developed to model proton 
conduction along one-dimensional water wires. Each lattice site represents a water 
molecule that can be in only one of three states; protonated, left-pointing, and 
right-pointing. Only a right (lcft)-pointing water can accept a proton from its 
left(right). Results of asymptotic mean field analysis and Monte-Carlo simulations 
for the three-species, open boundary exclusion model are presented and compared. 
The mean field results for the steady-state proton current suggest a number of 
regimes analogous to the low and maximal current phases found in the single 
species ASEP [B. Derrida, Physics Reports, 301, 65-83, (1998)]. We find that the 
mean field results arc accurate (compared with lattice Monte-Carlo simulations) 
only in the certain regimes. Refinements and extensions including more elaborate 
forces and pore defects are also discussed. 
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1. Introduction 

One-dimensional nonequilibrium processes have been used as models for many systems 
including interface growth [Q, || , traffic flow || || , pore transport M, [|, |7) , and quantum 
spin chains j8j. The simple asymmetric exclusion process (ASEP) has also been used 
to model mRNA translation |)|, and recent modifications have been applied in 
order to describe molecular motors and parallel pathways [jll], |l2|, E3. 

Building upon the exact and asymptotic results of Derrida, Evans, Rittenberg, 
Sandow, and others ||, [l4|, [l5|, |l7], |l8), the quantitative behavior of other more 
complex one-dimensional systems can be investigated. Asymmetric exclusion processes 
of multiple species have also been modelled, particularly for periodic boundary 
condition s [19| , E9], or in the context of spontaneous symmetry breaking for open 
systems [^1 . In the proton conduction problem, only one species (protons) arc 
transported while the other classes of particles correspond to different orientation 
states of the relatively fixed water molecules. 

In this paper, we present a restricted three-species exclusion process that models 
proton conduction through a "water wire." Our model is motivated by experiments on 
proton conduction through stable or transient water filaments suggesting conduction 
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Figure 1. (a) Schematic of an N = 11 three-species exclusion model that captures 
the steps in a Grotthuss mechanism of proton transport along a water wire. For 
typical gramicidin pores that span lipid membranes TV ~ 20 — 22. The transition 
rates are labeled in (b). Water dipolc kinks are denoted by thick lines. 



occurs via proton exchange between properly aligned water molecules 
"Grotthuss p4| " mechanism of proton conduction along a relatively static water 
chain has been proposed to give rise to the fast proton transfer rates relative to the 
permeation of distinguishable ions across structurally similar ion channels. 

2. Lattice Model 

Figure la is a schematic of model in which each "site" along the pore is occupied 
by an oxygen atom. This oxygen is either part of neutral water, H2O, or hydronium 
H3 + . We assume hydronium ions are spherically symmetric, in the sense that their 
three protons are either in a planar, left-right symmetric hybridization, or are rapidly 
smeared over the molecule. These singly protonated species are denoted "0" type 
particles. In contrast, the neutral waters have "permanent" dipole moments or electron 
lone-pair orientations that can nonetheless fluctuate. For simplicity, all water dipolcs 
(hydrogens) that point to right are classified as "+" particles, while those pointing 
left are "— " particles. Thus, each site can exist in only one of three states: 0, +, or 
— , corresponding to protonated, right, or left states, respectively. This simple choice 
of for labeling the occupancy configurations, <x; = {— 1,0, +1}, allows for fast integer 
computation in simulations. 

The transition rules are constrained by the orientation of the waters at each site 
and are defined in Fig. lb. A proton will enter the first site (i = 1) from the left 
reservoir and protonate the first water molecule with rate a only if the hydrogens 
of the first water is pointing to the right (such that its lone-pair electrons are left- 
pointing, ready to accept proton from the left reservoir). Similarly, if a proton exits 
from the first site into the left reservoir (with rate 7), it leaves the remaining hydrogens 
right-pointing. The analogous rules hold for the last site i = N with exit and entrance 
rates f3 and S respectively. The entrance rates a and 5 are functions to the proton 
concentration in the respective reservoirs as well as local electric potential gradients 
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near the pore mouths. In the pore interiors, a proton at site i can hop to the right(left) 
with rate p + (p_) only if the adjacent particle is a right(left)-pointing, unprotonated 
water molecule. If such a transition is made, the remaining water molecule at site i will 
be left(right) pointing. Physically, as a proton moves to the right, it leaves a decaying 
trail of — particles to its left. A left moving proton leaves a trail of + particles to its 
right. These trails of — or + particles are unable to accept another proton from the 
same direction. Clearly, protons can successively follow each other only if the waters 
can reorient without being protonated. The water reorientation rates are denoted k± 
(cf. Fig. lb). The rate- limiting steps in steady-state proton transfer across biological 
water channels is thought to be associated with water flipping. 

This water wire conduction model resembles the two species ASEPs treated by 



Evans et al. 21 if one identifies protons with vacancies, and allows the + and — 
elkass particles to interconvert with rates k±. Moreover, + particles here can only 
move left, and once they do, get converted to — type particles. Similarly, a — particle 
that moves to the right becomes a + particle. The allowed moves arc —0 <-> 0+. 
Although the master equation can be formally defined in terms of the probability 
Pjv(cti, tJ2, . . . ojv-i, cjy) for a given configuration {ct^}, I have failed to discover an 
analytical result for the steady state current proton (a — 0) current using matrix 
product methods ||, |l7], |l8). Since it is not always possible to find a finite 
dimensional matrix product representation, particularly for open boundaries, many 
practical applications, such as multipecies pore transport Q, must be treated using 
mean field approximations or numerical simulations. For the single species ASEP, the 
mean field result for the steady-state current is exact in the N — > oo limit, although the 
site occupancies do not agree with those from the exact solution. Therefore, a mean 
field theory in the present problem may shed light on the steady-states. We will obtain 
mean field analytic expressions that can be compared with numerical simulations. 

3. Mean field approximation 

For one-dimensional nearest lattice site hopping models, the mean-field assumption 
states that the probability that the system is in a configuration with occupancy <7j 
and cr;+i factors into a product of single particle probabilities: 

^2 Pn(ti,- ■ ■ ,(Ti,(T i+1 ,. . . ,T N ) = P 2 fa,a i+ l) ~ Pfa)Pfa+i), (1) 

n=o,±i 

where Pfa) = E CTj+1 =o,±i E{ Ti }=o,±i p n(ti, ■ . ■ , 0*, cr i+1 , . . . , t n ). The master 
equations thus become 



Pfa = 0) = -p+Pfa = 0)P(a i+1 = +1) -p-P(ai = 0)Pfa-i = -1) 

+ p+P(a i _ 1 = Q)P(a l = +l)+p_P{<j l+1 = Q)P(a q = -1) (2) 

Pfa = +1) = -k+Pia, = +1) + fe_P((7i = -1) 

-p+P(a ? = +l)P(o-i_i = 0) +p-P{<J l = 0)P(a i _ 1 = -1) (3) 

Pfa = -1) = k+Pfa = +1) - k-Pfa = -1) 

-p_Pfa +1 = 0)Pfa = -l)+p+Pfa = 0)Pfa +1 = +1), (4) 
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where P(oj — cr) are the individual probabilities that site j is occupied by a particle in 
state cr. Conservation at each site requires P(cr i = — 1) + P(<7j = 0) + P(<j i = +1) = 1. 
In the present treatment, we will assume that p±, k± are independent of neighboring 
configurations cr, i.e. that the transition rates are local and do not depend interaction 
with the states of neighboring sites. Within the mean field approximation, the 
instantaneous current of protons between sites i and i + 1 is 



J(i,i + 1) = p+P(<7< = 0)P{<r l+1 = +1) - p-P(a l+1 = 0)P{a t = -1) 
while the currents at the left and right boundaries are 
J(L, 1) = aP(cr 1 = +1) - 7 P(cr 1 = 0) 

J(N, R) = f3P{o- N = 0) - 5P(a N = -1). 



(5) 



(6) 



Upon defining steady-state particle occupancies over their probabilities, (P(oi = 
0)) = n, {P(cri = +1)) = Si, and (P(<Tj = —1)) = U. The steady-state currents are 
thus 



and 



( J(i, i + 1)) = p+riSi+i - p_r i+1 ti 



(J (L, 1)) = asi - jn = {a + 7)si +7*1-7 
(J(N, R)) = (3r N -8t N = p- (3s N - (/? + 5)t N . 



(7) 



(8) 



In steady-state, all currents have to be identical, J = (J(i,i + 1)) = (J(L,1)) = 
(J(N,R)). Thus, the constant, internal steady-state current (M) reads 



J = (J(i,i + 1)} = p+riS l+1 -p-r t+ iU 
while the equation for (P(<7; = +1)) = reads 

= k + Si - k-U -p-Ti+iU +p + r l s l+1 



k+Si — k-ti + J. 



(9) 



(10) 



From (|10|), we see that k + si — k-ti = — J is constant along the chain. Upon imposing 
particle conservation, r.; + Sj = 



tj = Ksi 



J K(i - n) 



J 



K 



k 



+ 



k- K + 1 k + + k- k- 
The mean end-site occupations can be found by substituting ( pj] ) into (| 

a K + 1 + a/k. 



ri 



rN 



J. 



a + j(K + 1) a + j(K + 1) 
K5 K+l + S/k 



-J. 



/3(K + 1) + K6 P(K + 1)+K6' 
For the internal sites, substitution of (|l^) into (||) yields 

= Ainn+i + A 3 r i+ i - A 2 n - Ai, 

where 



(11) 
yielding 

(12) 
(13) 



A x = J, A 2 



J 



P+ 



A, 



P- 



k + + J 



and A 



4 



k-P+ - fc+P- 



(14) 



(15) 
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Wc can now write ( 14 ) as a recursion relation of the form 
A 2 r t - A 1 



Upon redefining 



A 4 n + A 3 ' 

A 2 + A 3 _ M 
A 4 b A4' 



we find the recursion relation for £: 

c(J) 



+1 



1 



6 



c(J) 



A 2 A 3 + A 1 A 4 
(A 2 + A 3 )* • 



(16) 



(17) 



(18) 



Except for the more complicated dependence of c(J) on J, the recursion relation 
( fTcf ) is identical to that used in analyzing the mean-field properties of the totally and 
partially asymmetric exclusion process for a single species [[u], The solution to 
the difference equation is Q 

•ei^i-e+e-K^-e*- 1 ) (19) 



fc-i 



k-2 



<-k- 2 \ 



where 



1 ± v/1 - 4c( J) 



(20) 



Similarly, the boundary occupations £1 , £tv can be readily found by substituting (|12|) 
and (|l^) into (|l7]) and solving for £. 

5. i. Symmetric limit 

First consider the totally symmetric case where — fc+P-| < o(l/N) and the 

"unbiased" current vanishes as l/N. One important realisation of this condition is 
fc + //c_ = K = 1 and p + = p_, corresponding to waters without left-right preference 
along the wire, and to the absence of pondcrmotive forces on the wire protons. 
We name this limit the symmetric limit because the net bias on internal protons 
vanish [[]. The proton current is driven entirely by a proton number density difference 
between the two reservoirs, i.e. an asymmetry in injection and extraction rates at the 
boundaries. Upon iterating Eqn (Hq) with A4 cx k-p + — k + p- = 0, 



A, 



N-l 



Ai \ 1 - (A 2 /A 3 ) 
A 2 



N-l 



1-(A 2 /A a ) ■ 

From the boundary conditions at i = 1, N, we obtain an implicit equation for J; 



(21) 



M 

A, 



N-l 



K + 1 + a/k. 



J 



KS 



Ax 1 - (A 2 /Ai 



A, 



\N-1 



1 - (A 2 /A 3 ) 
K + l + S/k_ 
(3{K + 1) + K8 ^ P(K + l) + K5 



(22) 



J. 



f In the proton conduction literature, the term "symmetric" is often used to denote symmetric 
solutions, or equal proton concentrations in the two reservoirs. The current in this case is driven by 
an externally imposed transmembrane potential 
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where the Ai are functions of J given by Eqns. (|l5|). A closed analytic expression can 
be found as an asymptotic expansion about large N using the ansatz 



O-l CL2 / 1 



(23) 



where the ai are TV-independent coefficients that are functions of the kinetic rates. To 
lowest order, 



A 3 



N-l 



1 — d\ 



1 



Nk, 



N-l 



exp 



ax{K+l) 



At 0(1 /N°), (^) determines a\ and the steady-state current 

(3(K +1)+K5 7 (fc+ + + + fc_) 



J = 



fc + fc_ 



A(fc 4 



In 



7 (if + 1) + a /3(fc+ + fc_) + if<5(p_ + 



(24) 



+ 0(N- 2 ). (25) 



When /? = 7, the current (|25|) in the limit of small proton concentration differences 
can be expanded in powers of k-u — k + S, 



J 



/3p+(fc_a — k + S) 



+ O ({k^a - k+5) 2 ) + 0(1/N 2 ), 



(26) 



N [p(K + f ) + K5] (13 + S)(k + + k 

reflecting a current driven by the effective difference in boundary injection rates 
Finally, in the large a and 5 = limit, 

p_ \ jk + k-p 

17 



j 



(k+ + k-)N 



log 1 



aN(k+ + k-)(k+ +p-) 



+ 0(a- 2 N- 1 ). 



(27) 



3.2. Asymmetric limits 

For driven systems, where |fc_p+ — > 0(1/N), a finite current persists in 

the N — > oo limit. Mean field currents are found from the fixed points of the 
recursion relation (|lj). The average particle occupations at each site i can be then 
be reconstructed from the solution to £j. In analogy to the single-species partially 
asymmetric process, the condition c(Jm) = 1/4 leads to a single fixed point that 
determines Jm, via a quadratic equation which is solved by 



J 



2(p + k^ — p-k^ 



(P 



+ +P-) 



(P- +P+) 



+ fc_ + k + ± y/k + + k_ y/k- +k + +p + + 
when fi(J±) > i, ^(Jm) < |- 



P- 



(28) 



For the sake of simplicity, we assume fc_p + — fc + p_ > and furthermore {a, ft, 7, £, k±} 
are such that J > 0. Clearly J < (since J must be less than the smallest rate in 
the entire process) and A4 > 0. In the limit p + /(k + + fc_),p_/(fc + + k-) — > 0, only 
provides a physical solution, 



M 



P+k- ~P-*M 
k-L + fc_ 



P- 



8(k, 



k- 



+ 




(29) 



and the maximal current J M in this model 



For a purely asymmetric process p_ 
approaches the analogous maximal current expression of the single species ASEP, 




Figure 2. (a). Mean field phase diagram in the totally asymmetric limit 
5 = 7= p_ = 0. Three current regimes for p+ / (fc_|_ + fc_ ) = 10 are delineated 
by thick solid curves at values of a*/(fc+ + = 2/3*/(k+ + fc_) = \/l + x — 1. 
Thin dashed lines define the mean field phase diagram for fc+ = k— = 1.0. (b). 
The dependence of a*/(fc+ + on x = p+/(k+ + fc_). 



except for modification by the factor (K + representing the fraction of time the 
site to the immediate right of a proton is in the + state. 

Another possibility is c(J < J^) < 1/4 which leads to two fixed points, one 
stable at £_ and one unstable at £+. Following [JHJ and p8[ , we set in one case 
the left occupation £i = £_ to the stable value which sclf-consistently determines 
the current via c(J) = £i(l — £i), while equation ( |l2"| ) determines Thus 
c(J L ) = £i(J L )(l - Zi(Jl)) is solved by 

± = G(a,i,k±,p±) ^/G 2 (a,j,k±,p±) - iH(a,'y,k ± ,p ± )F(a,j,k ± ,p±) 
" L ~ 2F(a,j,k±,p±) 2F(a,~/,k±,p±) ' [ ' 

where 

F(x,y,k±,p±) = (k+ + k-)(x + k+ + k-)[y(p + +p-)+xp--U]U 



G(x,y,k±,p±) = {k+ + k-)U x 2 k-(k- - (xfc- - y(k+ + k-))U+ 

xy(2k + k- + 2k + p + + k-P- - k + p_ + 2k 2 _) + y 2 (k + + k-) 2 
H(x,y,k±,p±) = xyk-(k- + p+)U 2 (k±,p±) 



(32) 



U(k ± ,p ± ) = k_p + - k + p_. 

The injection rate-limited current ( |3l| ) depends only upon a, 7, k,p±. Upon iterating, 
we find the regime of validity for J ~ Jl is determined by £n < £+■ Thus, the currents 
can arise only when £i(J^ ) < 1/2 and ^,n{Jl) < 1 — Ci(^l) ano - correspond to a 
low proton density (r,) phase. 

Another limit of steady-state currents corresponds to the high density phase where 
£jv = £+, which is found from £jv(</i?)(l — £,n{Jr)) = c(Jr < Jjw) with 



r± = G(£,/3,fc T ,p T ) v /G 2 ((5,/3,fc T ,p T )-4g( ( 5,/3,fc T ,p T )F(^,/3,fc T ,p^) 
fi 2F(<5,/?,fc T ,p T ) 2F{5,p,k T ,p T ) 
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Equation (33) is also easily deduced from the solutions and the symmetry 



properties of c(J),r 1; and rjy under the transformation {a, 7, k±,p±, J} — > 
{<5,/3,A; T ,p T ,-J}. The solutions j| hold only if £ w (j|) > 1/2, &v(j|) > l-£i(jf ) 
and correspond to a high proton density phase. Although the currents represented 
by equations (28), (31), ([33]), and their associated regimes of validity define a 
potentially rich mean field current phase diagram, significant simplification, where 
simpler expressions can be obtained, is achieved by restricting ourselves to the totally 
asymmetric regime where 7 = S = p- = 0. In this case, A3 = 0, A4 = p+/(K + 1), 
and £ = fc_r/(fc_ — J). Upon using r\ and rjy given by ( [l2|) and (|l3|), the currents 
(^H), and (Js3"l), and their regimes of validity become 



J 



Af 



T:_ 



1 



2(fc+ + fc_) 



1 



1 



P+ 



k+ + k- 



ak— (1 — a/p+) 



where 



fc_(l-/3/ P+ )-/3fc + / P+ 

a* = (K + l)/3* = (fc+ + fc_) 



a < a*,/3 > — 



(A' + l) 
/9 < /3*,a > f3(K + l), 



(34) 



1 



k+ + k- 



1 



(35) 



The phase diagram defined by ([34j) and (pq) i s shown in Fig. 2a for x = p+/(fc + +fc_ ) = 
1, 10. In the k± — > 00 limit, the currents (|34|) reduce to modified forms of the standard 
currents in the single particle ASEP [pl| |l7f . 



4. Monte-Carlo simulations 



We now compare mean field results from both the symmetric and asymmetric 
limits with their corresponding steady-state currents obtained from lattice Monte 
Carlo (MC) simulations. The simulation was implemented by defining TV + 2 sites 
corresponding to the interior pore sites, plus the two reservoirs. At each time step, 
one of the TV + 2 sites was randomly chosen and a transition attempted. Transitions 
were accepted according to the Metropolis algorithm. At the next time step, another 
site was randomly chosen. The occupation at each reservoir end site was set to 
Co = o"at+i = so that protons were always ready to be injected with probability 
corresponding to the rates a, S. To obtain low noise steady-state values for the 
occupations and currents, instantaneous currents across all boundaries were averaged 
over O(10 8 — 10 9 ) time steps. 

Figures 3a, b show slices (for N = 50, k± = 0.1,1, and fixed a) through the 
phase diagram of the totally asymmetric case (Fig. 2). The mean field currents are 
represented by the thin solid and thick dashed curves for different values of a in each 
plot. Mean field theory seems to predict a first order transdition between the J7 and 
phases and a higher order transition to the phase. The curves labelled with 
O symbols represent steady-state currents computed from the MC simulations for the 
same fixed values of a. The mean field solutions always overestimate the currents 
with the error especially bad for small k± (compare Fig. 3a with Fig. 3b). This 
suggests correlations are suppressed when water flipping is fast, rendering the mean 
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Figure 3. Comparison of mean field and N = 50 Monte-Carlo simulation (filled 
circles) results for the steady-state current with 7 = 8 = p_ = and p+ = 2. (a) 
Currents as a function of /?/(fc_|- for fc_|_ = fc_ = 0.1 and a = 0.2,0.6. (b). 

Currents for fc+ = k— = 1 and a = 0.1,3. Thin solid curves correspond to the 
smaller values of a, while thick dashed curves depict currents with the larger a. 



field approximation more accurate. A similar qualitative trend arises in the symmetric 
case. 

Figure 4 shows the steady-state currents for N = 100 in the special symmetric 
case (p+ = p_ = 1) as a function of the equal flip rates k = k + = k-. The asymptotic 
mean-field currents for (a = 0.2,6 = 0.0) and (a = 1.0, 6 = 0.0) are plotted with 
thin and thick dashed curves, respectively. The mean-field results again overestimate 
the steady-state current J, and become increasingly accurate for large flip rates. The 
k — *• 00 limit in which mean field theory is asymptotically exact gives 



N 



a/3 -7<5 p- a 2 (3(/3 + 5) - ^7(7 + a) t I pt 



(2/3 + 5)(2 7 + a) k (2/3 + ,5)2(27 + a) 2 '~V*VJ" ^ 

Site occupations depicted in Fig. 5a clearly indicate the absence of boundary 
layers in the symmetric limit. The proton concentration is nearly linear along the 
chain and is fairly insensitive to the flipping rate k. Only for small k = 0.1 (thin lines, 
small symbols) does the profile deviate slightly from being linear. Not surprisingly, 
the mean field result for the steady-state current is substantially more accurate when 
k is large. Figure 5b shows the correlations 

S(vu <T i+ i) = ((1 - *m - of +1 )> - (1 - <jf)(l - of +1 ) (37) 

as functions of the lattice site i. We plot S(+, +) and ,9(0, +) (circles) for both k = 0.1 
(thin dashed lines) and k = 10(thick solid lines). Correlations S containing oi = — 
do not vary appreciably when k is changed. Note that £(+,+) for k = 0.1 is larger 
than 5(+,+) for k = 10. The larger flipping rates wash out the correlations and 



Proton 



10 



0.3- 



N = 100; p = y = 0.5, 8 = 0; p + = p_ = 1 




I I 

012345678 

k =k 
+ - 



Figure 4. Mean field (dashed) and simulated (solid) steady-state currents for a 
"symmetric" (fc— p+ = fc+p_) chain with N = 100, fi = 0.5,5 = 0,p+ = p_ = 1 
and a = 0.2, 1.0. 



suppresses S{+, +). Conversely, 5(0,+) for small k = 0.1 is negative since OH ► — 

is a relatively fast step that destroys the 0+ configuration. When k is large, this 
negative correlation is also supressed. 

Figures 6 contrasts concentration difference driven currents and potential 
difference driven currents. In figure 6a, 100 x J is plotted as a function of a for 
N = 20 and 6 = 0. The current is always sub-linear with respect to a and increases 
for small, increasing fc. For larger k, the current saturates as the rate limiting step 
becomes the forward hopping rate p. 

The site occupancies, averaged over the entire chain, are plotted in Fig. 6b. 
The curves plot a = A^EiliC 1 - °f)> ^- = (2A0 _1 EjLi - 1), and 

(j + = (2iV) _1 EiLi cr i( <7 i + 1); representing the chain-averaged proton, — type, and + 
type particle densities, respectively. As the proton entrance rate a is increased, the 
proton concentration in the pore increases at the expense of both — and + particles. 
However, since the current to the right increases, the average population of + particles 
decreases faster than that of the — particles. For larger k, the difference between — 
and + occupations vanishes as the high flip rate equalizes their probabilities; the only 
small remaining difference arises from the occasional passing of a proton converting 
+ to -. 

For potential driven flow, the forward and backward hopping rates were modified 
according to 

(a, ftp)-* (a,f3,p)e y / 2 

(38) 

(-y,S,q)^( 1 ,5,q)e- v / 2 

As a function of V, the current J can appear to be sublinear or super linear depending 
on the voltage regime (Fig. 6c). For small k and sufficiently large V, the rate limiting 
steps are the + <-» — flips and the current J is sublinear. The voltage-driven average 
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Figure 5. (a) Steady-state occupations of 0, — , + along the chain for a = 1,/? = 
7 = 0.5, 5 = and p± = 1. The larger symbols correspond to k± = k = 10, 
while the small, thin symbols indicate occupations associated with k = 0.1. 
The differences between en = and cr; = — are barely discernible, (b). The 
correlations +) and S(0, +) for small(fc = 0.1) and large(fc = 10) flip rates. 



particle densities for k = 0.1 (Fig. 6d) also show increasing proton concentrations and 
decreasing + particle densities. 



5. Conclusions 

We have formulated an interacting, three-species transport model (epitomized in 
Fig. 1) which is motivated by the Grotthuss mechanism [^4| of proton transfer 
along water wires. The dynamics resemble those of two-species ASEPs but with 
important differences that give rise to very different behavior. Our main results arc 
summarized by analytic expressions for the steady-state currents derived within a 



mean-field approximation. Currents in both symmetric (equations 25 and |2^) and 
asymmetric regimes (equations and [j3|) were found, with phases in the totally 

asymmetric limit (equation |34|) analogous to those of the standard ASEP. Monte-Carlo 
simulation agree with mean field analyses only in the limit of large flipping rates k± . 
Although qualitatively correct, the totally asymmetric mean field currents described 
by the phase diagram in Fig. 2, always overestimates the actual currents, as shown 
in Fig. 3. As apparent from Figs. 3,4, steady-state currents in both symmetric and 
asymmetric limits are very sensitive to subtle variations in average occupation and 
minute differences in nearest neighbor correlations plotted in Fig. 5. 

Our model does not include many important interactions that may be important 
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Figure 6. (a) J as a function of a for N = 20 and 6 = 0,p± = 1, and various k. 
The corresponding chain-averaged populations f or k = 0.1 are shown in (b). (c) 
Currents as a function of driving V defined by (J3q) . All J — V curves saturate 
(become sublinear and level off) at large V when flipping rates k become rate- 
limiting, (d) Averaged populations for k = 0.1 as functions of driving V. 



in real physical systems. Examples include dipole-dipolc interactions between adjacent 
waters, hydronium-hydronium repulsion between two adjacent protons, and potential 
induced ordering fields that preferentially orient the water dipoles. These effects can 
be readily incorporated by considering next-nearest neighbor interactions. The more 
complicated interaction energy can be used in weighting Monte-Carlo acceptance steps. 

A salient experimental finding in membrane potential-driven water wire proton 
conduction is the crossover from sublinear J — V curves to super linear curves as the pH 
in the two equal reservoirs is decreased |2^] . Multiple protons on the water wire have 
been implicated in the crossover to superlinear behavior. Unlike previous theories pa], 
our model allows for multiple proton entry into the water wire. In addition to finding 
a wider regime of supcrlincarity for large k±, we also found (not shown) stronger 
supcrlincarity when the injection rates (a, S) were large relative to extraction rates 
(P, 7). These results of the simple lattice model display a number of different regimes, 
provide relative estimates for fundamental tranition rates, and will serve as a starting 
point for more complicated models of measured J — V characteristics [^6| . 
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